#include <stdio.h>
#define emax 30
# define emaxp1 emax+1
#define dim 2*emax+1
#define gmax 1+emax/2

int hptest(int rt);
void makeword(int nq, int np);
void spairs(int e, int g);
void minpos(int np);
int hfirstmix(int nq, int np);
int hnextmix(int nq, int np);
int hplace(int i, int nq, int np);
int hmove(int i);
int colour(int np);
int htest(int rt);
int hpcount(int rt);
int hcount(int rt);
void count(int rt);
void pcount(int rt);
void test(int rt);
void ptest(int rt);
int Lehman(int e, int v, int vert[dim], int next[dim], int r, int code[dim], int c);
int Alfred(int e, int v, int vert[dim], int next[dim], int r, int code[dim], int c);
void decode(int e, int v, int vert[dim], int next[dim]);
int nextsys(int np, int nq);
void firstsys(int np, int nq);
int ppairs(int n);
void qpairs(int n);
int genus(int n);
void firstint(int n); int nextint(int n);
void firstpar(int n);
int nextpar(int n);
void firstbrac(int n);
int nextbrac(int n);
int bpairs(int n);

int p[dim], q[dim], dir[emaxp1], qpos1[emaxp1], qpos2[emaxp1], col[emaxp1], wall[dim],
	qmate[dim], pmate[dim], smate[dim], qpos[dim], s[dim], pr[dim], prpos[dim], bug;

int main(void)
{
	int tst, rt, allg, hyper;
	
	while (1)
	{
		printf("Enter 1 to test the program, 2 to count maps or hypermaps, 0 to quit: ");
		scanf("%d",&tst);
		if (tst==0) break;
		printf("Enter 1 for hypermaps, 0 for maps: ");
		scanf("%d",&hyper);
		if (hyper>0)
		{
			printf("Enter 0 for rooted hypermaps, 1 for sensed hypermaps, 2 for unsensed hypermaps: ");
			scanf("%d",&rt);
			printf("Enter 0 for planar hypermaps or 1 for hypermaps of all genus: ");
			scanf("%d",&allg);
			if (tst==1) { if (allg) htest(rt); else hptest(rt); }
			else { if (allg) hcount(rt); else hpcount(rt); }
		}
		else
		{
			printf("Enter 0 for rooted maps, 1 for sensed maps, 2 for unsensed maps: ");
			scanf("%d",&rt);
			printf("Enter 0 for planar maps or 1 for maps of all genus: ");
			scanf("%d",&allg);			
			if (tst==1) { if (allg) test(rt); else ptest(rt); }
			else { if (allg) count(rt); else pcount(rt); }
		}
	}
	printf("End of execution.\n");
	return 0;
}

int hptest(int rt)
{
	int i, j, np, nq, flag, nw, unroot, d, r,max,vert[dim],next[dim],code[dim],prev[dim];
	char cp;
	
	if (rt==0) printf("Testing canonical mixing of pairs of Dyck words.\n");
	else printf("Testing Lehman coding and decoding and maximality.\n");
	if (rt>1) printf("Also testing maximality compared to reflected map.\n");
	while (1)
	{
		printf("Enter number of parenthesis pairs <=%d, -1 to quit: ",emax);
		scanf("%d",&np);
		if (np<0) break;
		printf("Enter a balanced parenthesis system with %d pair",np);
		if (np!=1) printf("s");
		printf(".\n");
		flag=0;
		for (i=1;i<=2*np;i++)
		{
			do scanf("%c",&cp); while (cp!='('&&cp!=')');
			if (cp=='(') p[i]=0;
			else if (cp==')') p[i]=-1;
			else flag=1;
		}
		p[0]=-1;
		if (!flag) flag=ppairs(np);
		if (flag) {printf("Invalid parenthesis system.\n"); continue;}
		printf("    ");
		for (i=1; i<=2*np; i++)
		{
			if (p[i]==0) printf("  ("); else printf("  )");
		}
		printf("\n ");
		for (i=0; i<=2*np; i++) printf("%3d",i);
		printf("\nmate");
		for (i=1; i<=2*np; i++) printf("%3d",pmate[i]);
		printf("\nw");
		minpos(np);
		for (i=0; i<=2*np; i++) printf("%3d",wall[i]);
		nw=colour(np);
		printf("\n");
		printf("Vertex  ");
		for (i=1; i<=np+1; i++) printf("%3d",i);
		printf("\n");
		printf("Colour  ");
		for (i=1; i<=np+1; i++)
		{
			if (col[i]==1) printf("  w");
			else printf("  b");
		}
		printf("\n");
		printf("Number of white vertices = %d, number of black vertices = %d\n",nw,np+1-nw);
		printf("Enter number of bracket pairs <= %d, -1 to quit: ",emax-np);
		scanf("%d",&nq);
		if (nq<0) break;
		if (nq>emax-np) continue;
		printf("Enter a balanced bracket system with %d pair",nq);
		if (nq!=1) printf("s");
		printf(".\n");
		flag=0;
		for (i=0;i<2*nq;i++)
		{
			do scanf("%c",&cp); while (cp!='['&&cp!=']');
			if (cp=='[') q[i]=2;
			else if (cp==']') q[i]=1;
			else flag=1;
		}
		q[2*nq]=0;
		if (!flag) flag=bpairs(nq);
		if (flag) {printf("Invalid bracket system.\n"); continue; }
		for (i=0; i<2*nq; i++)
		{
			if (q[i]==2 )printf("  ["); else printf("  ]");
		}
		printf("\n");
		for (i=0; i<2*nq; i++) printf("%3d",i);
		printf("\n");
		for (i=0; i<2*nq; i++) printf("%3d",qmate[i]);
		printf("\n");
		printf("Here are all the ways to mix those two words.\n");
		flag=hfirstmix(nq,np); /* Mixing these two words in all possible ways. */
		while (!flag)
		{
			for (j=0;j<2*nq;j++) printf("%3d",qpos[j]);
			makeword(nq,np);
			for (j=1;j<=2*(nq+np);j++)
			{
				if (s[j]==2) printf("  [");
				else if (s[j]==1) printf("  ]");
				else if (s[j]==0) printf("  (");
				else if (s[j]==-1) printf("  )");
				else printf(" Oops!" );
			}
			printf("\n");
			if (rt>0)
			{
				printf("Enter 1 to decode this system or 0 not to: ");
				scanf("%d",&unroot);
				if (unroot)
				{
					spairs(np+nq, 0);
					for (i=1;i<=2*(np+nq);i++) printf("%3d",i);
					printf("\n");
					for (i=1;i<=2*(np+nq);i++) printf("%3d",smate[i]);
					printf("\n");
					decode(np+nq, np+1, vert, next);
					printf("vert: ");
					for (i=1; i<=2*(np+nq); i++) printf("%3d",vert[i]);
					printf("\n");
					printf("next: ");
					for (i=1; i<=2*(np+nq); i++) printf("%3d",next[i]);
					printf("\n");
					i = Alfred(np+nq, np+1, vert, next, 1, code, 0);
					printf("And here is the code of the rooted map.\n");
					for (i=1;i<=2*(np+nq);i++)
					{
						printf("%3d",code[i]);
						if (code[i]!=s[i])
						{
							printf("  Oops!\n"); break;
						}
					}
					printf("\n");
					printf("Enter 1 to test this code for maximality or 0 not to: ");
					scanf("%d",&unroot);
					if (unroot)
					{
						max=1;
						for (r=2; r<=2*(np+nq); r++) if (col[vert[r]]==1)
						{
							d = Alfred(np+nq, np+1, vert, next, r, code, 1);
							if (d>0)
							{
								printf("This code is not maximal.  Rooting at %d is bigger.\n",r);
								i = Alfred(np+nq, np+1, vert, next, r, code, 0);
								for (i=1; i<=2*(np+nq); i++) printf("%3d", code[i]);
								printf("\n");
								max=0;
								break;
							}
						}
						if (max)
						{
							printf("This code is maximal.\n");
							if (rt>1)
							{
								for (i=1;i<=2*(np+nq);i++) prev[next[i]]=i;
								printf("vert: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",vert[i]);
								printf("\n");
								printf("next: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",next[i]);
								printf("\n");
								printf("prev: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",prev[i]);
								printf("\n");
								for (r=1; r<=2*(np+nq); r++) if (col[vert[r]]==1)
								{
									d = Alfred(np+nq, np+1, vert, prev, r, code, 1);
									if (d>0)
									{
										printf("This code is not reflect-maximal.  Rooting at %d is bigger.\n",r);
										i = Alfred(np+nq, np+1, vert, prev, r, code, 0);
										for (i=1; i<=2*(np+nq); i++) printf("%3d", code[i]);
										printf("\n");
										max=0;
										break;
									}
								}
								if (max) printf("This code is reflect-maximal.\n");
							} /* end if reflect-maximality test flag */
						} /* end if max */
					} /* end if maximality test flag */
				} /* end if coding flag */
			} /* end if rt>0 */
			flag=hnextmix(nq,np);
		} /* end while (!flag) */
		printf("These are the only ways to mix those two words.\n");
	} /* end while(1) */
	return 0;
}

void makeword(int nq, int np)
/* makes a shuffle of parentheses and brackets or integers */
{
	int i, ip, iq;
	
	if (nq<=0) for (i=1;i<=2*np;i++) s[i]=p[i];
	else
	{
		ip=1;
		iq=0;
		for (i=1;i<=2*(nq+np);i++)
		{
			if (iq<2*nq&&qpos[iq]<ip)
			{
				s[i]=q[iq];
				iq++;
			}
			else
			{
				s[i]=p[ip];
				ip++;
			}
		} /* end for */
	}
} /* end makeword */

void spairs(int e, int g)
/* Finds the position smate[i] of the mate of s[i]  */
{
	int i, stk1[dim], stk2[dim], ind1, ind2;
	
	stk1[0]=0;
	stk2[0]=0;
	ind1=0;
	ind2=0;
	if (g>0) for (i=1;i<=e;i++) stk2[i]=0;
	for (i=1;i<=2*e;i++)
	{
		if (s[i]==0) /* left parenthesis; push */
		{
			ind1++;
			stk1[ind1]=i;
		}
		else if (s[i]==-1) /* right parenthesis; pop */
		{
			smate[i]=stk1[ind1];
			smate[stk1[ind1]]=i;
			ind1--;
		}
		else /* bracket or integer */
		{
			if (g==0) /* parenthesis-bracket system */
			{
				if (s[i]==2) /* left bracket; push */
				{
					ind2++;
					stk2[ind2]=i;
				}
				else /* right bracket; pop */
				{
					smate[i]=stk2[ind2];
					smate[stk2[ind2]]=i;
					ind2--;
				}
			}	/* end parenthesis-bracket */
			else /* parenthesis-integer system */
			{
				if (stk2[s[i]]==0) stk2[s[i]]=i; /* first occurrence */
				else							/* second occurrence */
				{
					smate[i]=stk2[s[i]];
					smate[stk2[s[i]]]=i;
				}
			} /* end parenthesis-integer */
		} /* end outer else */
	} /* end for */
} /* end spairs */

int hfirstmix(int nq, int np)
{
	int i, flag;
	i=2*nq-1;
	if (nq==0) return 0;
	flag=hplace(i,nq,np);
	while (1)
	{
		if (flag)
		{
			i++;
			if (i>=2*nq) return 1;
			flag=hmove(i);
		}
		else
		{
			i--;
			if (i<0) return 0;
			flag=hplace(i,nq,np);
		}
	}
} /* end hfirstmix */

int hnextmix(int nq, int np)
{
	int i, flag;

	if (nq==0) return 1;
	i=0;
	flag=hmove(i);
	while (1)
	{
		if (flag)
		{
			i++;
			if (i>=2*nq) return 1;
			flag=hmove(i);
		}
		else
		{
			i--;
			if (i<0) return 0;
			flag=hplace(i,nq,np);
		}
	}		
} /* end hnextmix */	

int hplace(int i, int nq, int np)
{
	int pl;
	if (qmate[i]>i) /* first occurrence */
	{
		if ((qpos[qmate[i]]-qpos[i+1])%2==0) pl=qpos[i+1]-1;
		/* There must be an odd number of parentheses between matching brackets or integers. */
		else pl=qpos[i+1];
		if (pl<wall[qpos[qmate[i]]]) return 1; /* It would make the forbidden subword. 
*/
	}
	else
	{
		if (i==2*nq-1) pl=2*np; /* put the rightmost bracket or integer in the rightmost slot*/
		else pl=qpos[i+1];
		if (pl<2) return 1;
	}
	qpos[i]=pl;
	return 0;
} /* end hplace */

int hmove(int i)
{
	if (qmate[i]>i) /* first occurrence */
	{
		if (qpos[i]>wall[qpos[qmate[i]]]+1) {qpos[i]--;qpos[i]--;}
		else return 1;
	}
	else
	{
		if (qpos[i]>2) qpos[i]--;
		else return 1;
	}
	return 0;
} /* end hmove */

void minpos(int np)
{
	int ist, i, stk[emaxp1];
	stk[0]=0; ist=0; wall[2*np]=0;
	for (i=2*np; i>0; i--)
	{
		if (p[i]<0) /* right parenthesis */
		{
			ist++;
			stk[ist]=pmate[i]; /* push */
		}
		else ist--; /* pop */
		wall[i-1]=stk[ist];
	} /* end for */
} 

int colour(int np)
{
	int i, cl, iv, nw;
	
	col[1]=1; /* 1 = white, 0 = black */
	cl=1; /* colour of current vertex */
	iv=1; /* index of current new vertex */
	nw=1; /* number of white vertices */
	for (i=1; i<=2*np; i++)
	{
		if (cl==1) cl=0; else cl=1; /* change colour */
		if (p[i]==0) /* left parenthesis, new vertex */
		{
			iv++;
			col[iv]=cl;
			if (cl==1) nw++;
		} /* end if */
	} /* end for */
	return nw;
} /* end colour */

int htest(int rt)
{
	int i, j, np, nq, flag, nw, unroot, d, r,max,occ[emaxp1],vert[dim],next[dim],code[dim],prev[dim];
	char cp;
	
	if (rt==0) printf("Testing canonical mixing of Dyck word with integer system.\n");
	else printf("Testing Lehman coding and decoding and maximality.\n");
	if (rt>1) printf("Also testing maximality compared to reflected map.\n");
	while (1)
	{
		printf("Enter number of parenthesis pairs <=%d, -1 to quit: ",emax);
		scanf("%d",&np);
		if (np<0) break;
		printf("Enter a balanced parenthesis system with %d pair",np);
		if (np!=1) printf("s");
		printf(".\n");
		flag=0;
		for (i=1;i<=2*np;i++)
		{
			do scanf("%c",&cp); while (cp!='('&&cp!=')');
			if (cp=='(') p[i]=0;
			else if (cp==')') p[i]=-1;
			else flag=1;
		}
		p[0]=-1;
		if (!flag) flag=ppairs(np);
		if (flag) {printf("Invalid parenthesis system.\n"); continue;}
		printf("    ");
		for (i=1; i<=2*np; i++)
		{
			if (p[i]==0) printf("  ("); else printf("  )");
		}
		printf("\n ");
		for (i=0; i<=2*np; i++) printf("%3d",i);
		printf("\nmate");
		for (i=1; i<=2*np; i++) printf("%3d",pmate[i]);
		printf("\nw");
		minpos(np);
		for (i=0; i<=2*np; i++) printf("%3d",wall[i]);
		nw=colour(np);
		printf("\n");
		printf("Vertex  ");
		for (i=1; i<=np+1; i++) printf("%3d",i);
		printf("\n");
		printf("Colour  ");
		for (i=1; i<=np+1; i++)
		{
			if (col[i]==1) printf("  w");
			else printf("  b");
		}
		printf("\n");
		printf("Number of white vertices = %d, number of black vertices = %d\n",nw,np+1-nw);
		printf("Enter number of integer pairs <= %d, -1 to quit: ",emax);
		scanf("%d",&nq);
		if (nq<0) break;
		if (nq>emax) continue;
		printf("Enter each integer from 1 to %d twice with the first occurrences in increasing order.\n",nq);
		flag=0;
		for (i=1;i<=nq;i++) occ[i]=0;
		for (i=0;i<2*nq;i++)
		{
			scanf("%d",&q[i]);
			occ[q[i]]++;
			if (occ[q[i]]==1) qpos1[q[i]]=i; else if (occ[q[i]]==2) qpos2[q[i]]=i; else flag=1;
		}
		q[2*nq]=0;
		for (i=1; i<=nq; i++)
			if (occ[i]!=2) {flag=1; break;}
		if (qpos1[1]!=0)  flag=1;
		for (i=2; i<=nq; i++)
			if (qpos1[i]<=qpos1[i-1]) {flag=1; break;}
		if (flag) {printf("Invalid integer system\n"); continue;}
		qpairs(nq);
		for (i=0; i<2*nq; i++) printf("%3d",q[i]);
		printf("   ");
		for (i=0; i<2*nq; i++) printf("%3d",qmate[i]);
		printf("  genus = %d\n", genus(nq));
		printf("Here are all the ways to mix those two words.\n");
		flag=hfirstmix(nq,np); /* Mixing these two words in all possible ways. */
		while (!flag)
		{
			for (j=0;j<2*nq;j++) printf("%3d",qpos[j]);
			makeword(nq,np);
			for (j=1;j<=2*(nq+np);j++)
			{
				if (s[j]==0) printf("  (");
				else if (s[j]==-1) printf("  )");
				else printf("%3d",s[j] );
			}
			printf("\n");
			if (rt>0)
			{
				printf("Enter 1 to decode this system or 0 not to: ");
				scanf("%d",&unroot);
				if (unroot)
				{
					spairs(np+nq, 1);
					for (i=1;i<=2*(np+nq);i++) printf("%3d",i);
					printf("\n");
					for (i=1;i<=2*(np+nq);i++) printf("%3d",smate[i]);
					printf("\n");
					decode(np+nq, np+1, vert, next);
					printf("vert: ");
					for (i=1; i<=2*(np+nq); i++) printf("%3d",vert[i]);
					printf("\n");
					printf("next: ");
					for (i=1; i<=2*(np+nq); i++) printf("%3d",next[i]);
					printf("\n");
					i = Lehman(np+nq, np+1, vert, next, 1, code, 0);
					printf("And here is the code of the rooted map.\n");
					for (i=1;i<=2*(np+nq);i++)
					{
						printf("%3d",code[i]);
						if (code[i]!=s[i])
						{
							printf("  Oops!\n"); break;
						}
					}
					printf("\n");
					printf("Enter 1 to test this code for maximality or 0 not to: ");
					scanf("%d",&unroot);
					if (unroot)
					{
						max=1;
						for (r=2; r<=2*(np+nq); r++) if (col[vert[r]]==1)
						{
							d = Lehman(np+nq, np+1, vert, next, r, code, 1);
							if (d>0)
							{
								printf("This code is not maximal.  Rooting at %d is bigger.\n",r);
								i = Lehman(np+nq, np+1, vert, next, r, code, 0);
								for (i=1; i<=2*(np+nq); i++) printf("%3d", code[i]);
								printf("\n");
								max=0;
								break;
							}
						}
						if (max)
						{
							printf("This code is maximal.\n");
							if (rt>1)
							{
								for (i=1;i<=2*(np+nq);i++) prev[next[i]]=i;
								printf("vert: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",vert[i]);
								printf("\n");
								printf("next: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",next[i]);
								printf("\n");
								printf("prev: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",prev[i]);
								printf("\n");
								for (r=1; r<=2*(np+nq); r++) if (col[vert[r]]==1)
								{
									d = Lehman(np+nq, np+1, vert, prev, r, code, 1);
									if (d>0)
									{
										printf("This code is not reflect-maximal.  Rooting at %d is bigger.\n",r);
										i = Lehman(np+nq, np+1, vert, prev, r, code, 0);
										for (i=1; i<=2*(np+nq); i++) printf("%3d", code[i]);
										printf("\n");
										max=0;
										break;
									}
								}
								if (max) printf("This code is reflect-maximal.\n");
							} /* end if reflect-maximality test flag */
						} /* end if max */
					} /* end if maximality test flag */
				} /* end if coding flag */
			} /* end if rt>0 */
			flag=hnextmix(nq,np);
		} /* end while (!flag) */
		printf("These are the only ways to mix those two words.\n");
	} /* end while(1) */
	return 0;
}

int hpcount(int rt)
{
	int i, v, e, np, nq, emx, emn, vmx, vmn, prt, r, vert[dim], next[dim], code[dim], prev[dim], max;
	int nw, nb, flag;
	long int ct, root[dim], sensed[dim][dim], unsensed[dim][dim], sum, sums, sumu;
	
	printf("Counting rooted");
	if (rt==1) printf(" and sensed");
	if (rt>=2) printf(", sensed and unsensed");
	printf(" planar hypermaps.\n");
	while (1)
	{
		printf("\nEnter maximum number of darts <= %d or -1 to quit: ",emax);
		scanf("%d",&emx);
		if (emx<0) break;
		if (emx>emax) continue;
		printf("Enter minimum number of darts between 1 and %d: ",emx);
		scanf("%d",&emn);
		if (emn>emx) emn=emx;
		if (emn<1) emn=1;
		for (e=emn; e<=emx; e++)
		{
			printf("\nNumber of darts = %d\n",e);
			printf("Enter maximum total number of edges and vertices between 2 and %d: ",e+1);
			scanf("%d",&vmx);
			if (vmx>e+1) vmx=e+1;
			if (vmx<2) vmx=2;
			printf("Enter minimum total number of edges and vertices between 2 and %d: ",vmx);
			scanf("%d",&vmn);
			if (vmn<2) vmn=2;
			printf("Enter 1 to print the Lehman code for each rooted hypermap or 0 not to: ");
			scanf("%d",&prt);
			printf("\nNumber of rooted planar hypermaps with d darts, v vertices and e edges.\n");
			printf("  d   v   e\n");
			sum=0;
			if (rt>0) sums=0;
			if (rt>1) sumu=0;
			for (v=vmn; v<=vmx; v++)
			{
				np=v-1;
				nq=e-np;
				for (i=0; i<=vmx; i++) root[i]=0;
				ct=0;
				if (rt>0) for (i=0; i<=vmx; i++) sensed[v][i]=0;
				if (rt>1) for (i=0; i<=vmx; i++) unsensed[v][i]=0;
				firstpar(np);
				do
				{
					if (ppairs(np)>0) {printf("Invalid parenthesis system.\n"); break;}
					minpos(np);
					nw=colour(np);
					nb=v-nw;
					firstbrac(nq);
					do
					{
						if (bpairs(nq)>0) {printf("Invalid bracket system.\n"); break;}

						flag=hfirstmix(nq,np); /* Mixing these two words in all possible ways. */
						while (!flag)
						{
							makeword(nq,np);
							ct++;
							root[nw]++;
							sum++;
							if (prt)
							{
								printf("nw=%3d, nb=%3d, ct=%5ld:",nw,nb,ct);
								for (i=1; i<=2*e; i++)
								{
									if (s[i]==0) printf(" (");
									else if (s[i]==-1) printf(" )");
									else if (s[i]==2) printf(" [");
									else if (s[i]==1) printf(" ]");
									else {printf(" Oops!\n"); return 0;}
								} /* end for i */
							} /* end if(prt) */
							if (rt>0)
							{
								spairs(e, 0);
								decode(e,v,vert,next);
								max=1;
								for (r=2;r<=2*e;r++) if (col[vert[r]])
								{
									i = Alfred(e, v, vert, next, r, code, 1);
									if (i>0)
									{
										max=0;
										break;
									}
								} /* end for r */
								if (max)
								{
									sensed[v][nw]++;
									sums++;
									if (prt) printf("  Sensed map #%ld.",sensed[v][nw]);
									if (rt>1)
									{
										for (i=1;i<=2*e;i++) prev[next[i]]=i;
										for (r=1;r<=2*e;r++) if (col[vert[r]])
										{
											i = Alfred(e, v, vert, prev, r, code, 1);
											if (i>0)
											{
												max=0;
												break;
											}
										} /* end for r */
										if (max)
										{
											unsensed[v][nw]++;
											sumu++;
											if (prt) printf(" Unsensed map #%ld.",unsensed[v][nw]);
										} /* end if reflex-maximality test */
									} /* end if rt>1 */
								} /* end if maximality test */
							} /* end if rt>0 */
							if (prt) printf("\n");
							flag=hnextmix(nq,np);
						}/* end while(!flag) */
					} 
					while (nextbrac(nq)>0);
				}
				while (nextpar(np)>0);
				for (nw=1;nw<v;nw++)
				{
					printf("%3d%4d%4d",e,nw,v-nw);
					printf("%12ld\n",root[nw]);
				}
			} /* end for v */
			printf("%3d     sum",e);
			printf("%12ld\n",sum);
			if (rt>0)
			{
				printf("\nNumber of sensed planar hypermaps with d darts, v vertices and e edges.\n");
				printf("  d   v   e\n");
				for (v=vmn; v<=vmx; v++)
				{
					for (nw=1;nw<v;nw++)
					{
						printf("%3d%4d%4d",e,nw,v-nw);
						printf("%12ld\n",sensed[v][nw]);
					}
				}
				printf("%3d     sum",e);
				printf("%12ld\n",sums);				
			} /* end if rt>0 */
			if (rt>1)
			{
				printf("\nNumber of unsensed planar hypermaps with d darts, v vertices and e edges.\n");
				printf("  d   v   e\n");
				for (v=vmn; v<=vmx; v++)
				{
					for (nw=1;nw<v;nw++)
					{
						printf("%3d%4d%4d",e,nw,v-nw);
						printf("%12ld\n",unsensed[v][nw]);
					}
				}
				printf("%3d     sum",e);
				printf("%12ld\n",sumu);				
			} /* end if rt>1 */			
		} /* end for e*/
	} /* end while */
	printf("End of execution of hpcount.\n");
	return 0;
} /* end hpcount */

int hcount(int rt)
{
	int i, g, v, e, np, nq, emx, emn, vmx, vmn, prt, r, vert[dim], next[dim], code[dim], prev[dim], max;
	int nw, nb, flag;
	long int ct, cts[dim][dim], ctu[dim][dim], root[gmax][dim], sensed[dim][gmax][dim], unsensed[dim][gmax][dim],
	sum[gmax], sumsum, sums[gmax], sumsums, sumu[gmax], sumsumu, ctr[dim];
	
	
	printf("Counting rooted");
	if (rt==1) printf(" and sensed");
	if (rt>=2) printf(", sensed and unsensed");
	printf(" hypermaps.\n");
	while (1)
	{
		printf("\nEnter maximum number of darts <= %d or -1 to quit: ",emax);
		scanf("%d",&emx);
		if (emx<0) break;
		if (emx>emax) continue;
		printf("Enter minimum number of darts between 1 and %d: ",emx);
		scanf("%d",&emn);
		if (emn>emx) emn=emx;
		if (emn<1) emn=1;
		for (e=emn; e<=emx; e++)
		{
			printf("\nNumber of darts = %d\n",e);
			printf("Enter maximum total number of edges and vertices between 2 and %d: ",e+1);
			scanf("%d",&vmx);
			if (vmx>e+1) vmx=e+1;
			if (vmx<2) vmx=2;
			printf("Enter minimum total number of edges and vertices between 2 and %d: ",vmx);
			scanf("%d",&vmn);
			if (vmn<2) vmn=2;
			printf("Enter 1 to print the Lehman code for each rooted hypermap or 0 not to: ");
			scanf("%d",&prt);
			printf("\nNumber of rooted hypermaps with d darts, v vertices and e edges.\n");
			printf("  d   v   e");
			for (g=0;g<=(e-1)/2;g++) printf("         g=%1d",g);
			printf("       all g\n");
			sumsum=0;
			for (g=0;g<=e/2;g++) sum[g]=0;
			if (rt>0)
			{
				sumsums=0;
				for (g=0;g<=e/2;g++) sums[g]=0;
			}
			if (rt>1)
			{
				sumsumu=0;
				for (g=0;g<=e/2;g++) sumu[g]=0;
			}
			for (v=vmn; v<=vmx; v++)
			{
				np=v-1;
				nq=e-np;
				ct=0;
				for (i=0; i<=vmx; i++)
				{
					ctr[i]=0;
					for (g=0; g<=vmx; g++) root[g][i]=0;
					if (rt>0) for (g=0; g<=vmx; g++)
					{
						sensed[v][g][i]=0;
						cts[v][i]=0;
					}
					if (rt>1) for (g=0; g<=vmx; g++)
					{
						unsensed[v][g][i]=0;
						ctu[v][i]=0;
					}
				}
				firstpar(np);
				do
				{
					if (ppairs(np)>0) {printf("Invalid parenthesis system.\n"); break;}
					minpos(np);
					nw=colour(np);
					nb=v-nw;
					firstint(nq);
					do
					{
						qpairs(nq);
						g=genus(nq);
						flag=hfirstmix(nq,np); /* Mixing these two words in all possible ways. */
						while (!flag)
						{
							makeword(nq,np);
							ct++;
							ctr[nw]++;
							root[g][nw]++;
							sum[g]++;
							sumsum++;
							if (prt)
							{
								printf("nw=%3d, nb=%3d, ct=%5ld:",nw,nb,ct);
								for (i=1; i<=2*e; i++)
								{
									if (s[i]==0) printf(" (");
									else if (s[i]==-1) printf(" )");
									else {printf("%3d",s[i]);}
								} /* end for i */
								printf("  genus = %d",g);
							} /* end if(prt) */
							if (rt>0)
							{
								spairs(e, 1);
								decode(e,v,vert,next);
								max=1;
								for (r=2;r<=2*e;r++) if (col[vert[r]])
								{
									i = Lehman(e, v, vert, next, r, code, 1);
									if (i>0)
									{
										max=0;
										break;
									}
								} /* end for r */
								if (max)
								{
									cts[v][nw]++;
									sensed[v][g][nw]++;
									sums[g]++;
									sumsums++;
									if (prt) printf("  Sensed map #%ld.",sensed[v][g][nw]);
									if (rt>1)
									{
										for (i=1;i<=2*e;i++) prev[next[i]]=i;
										for (r=1;r<=2*e;r++) if (col[vert[r]])
										{
											i = Lehman(e, v, vert, prev, r, code, 1);
											if (i>0)
											{
												max=0;
												break;
											}
										} /* end for r */
										if (max)
										{
											ctu[v][nw]++;
											unsensed[v][g][nw]++;
											sumu[g]++;
											sumsumu++;
											if (prt) printf(" Unsensed map #%ld.",unsensed[v][g][nw]);
										} /* end if reflex-maximality test */
									} /* end if rt>1 */
								} /* end if maximality test */
							} /* end if rt>0 */
							if (prt) printf("\n");
							flag=hnextmix(nq,np);
						}/* end while(!flag) */
					} 
					while (nextint(nq)>0);
				}
				while (nextpar(np)>0);
				for (nw=1;nw<v;nw++)
				{
					printf("%3d%4d%4d",e,nw,v-nw);
					for (g=0;g<=(e-1)/2;g++) printf("%12ld",root[g][nw]);
					printf("%12ld\n",ctr[nw]);
				}
			} /* end for v */
			printf("%3d     sum",e);
			for (g=0;g<=(e-1)/2;g++) printf("%12ld",sum[g]);
			printf("%12ld\n",sumsum);
			if (rt>0)
			{
				printf("\nNumber of sensed hypermaps with d darts, v vertices and e edges.\n");
				printf("  d   v   e");
				for (g=0;g<=(e-1)/2;g++) printf("         g=%1d",g);
				printf("        all g\n");
				for (v=vmn; v<=vmx; v++)
				{
					for (nw=1;nw<v;nw++)
					{
						printf("%3d%4d%4d",e,nw,v-nw);
						for (g=0;g<=(e-1)/2;g++) printf("%12ld",sensed[v][g][nw]);
						printf("%12ld\n",cts[v][nw]);
					}
				}
				printf("%3d     sum",e);
				for (g=0;g<=(e-1)/2;g++) printf("%12ld",sums[g]);
				printf("%12ld\n",sumsums);
			} /* end if rt>0 */
			if (rt>1)
			{
				printf("\nNumber of unsensed hypermaps with d darts, v vertices and e edges.\n");
				printf("  d   v   e");
				for (g=0;g<=(e-1)/2;g++) printf("         g=%1d",g);
				printf("        all g\n");
				for (v=vmn; v<=vmx; v++)
				{
					for (nw=1;nw<v;nw++)
					{
						printf("%3d%4d%4d",e,nw,v-nw);
						for (g=0;g<=(e-1)/2;g++) printf("%12ld",unsensed[v][g][nw]);
						printf("%12ld\n",ctu[v][nw]);
					}
				}
				printf("%3d     sum",e);
				for (g=0;g<=(e-1)/2;g++) printf("%12ld",sumu[g]);
				printf("%12ld\n",sumsumu);
			} /* end if rt>1 */			
		} /* end for e*/
	} /* end while */
	printf("End of execution of hcount.\n");
	return 0;
}

void count(int rt)
{
	int i, g, v, e, np, nq, emx, emn, vmx, vmn, prt, r, vert[dim], next[dim], code[dim], prev[dim], max;
	long int ct, cts[dim], ctu[dim], root[gmax], sensed[dim][gmax], unsensed[dim][gmax],
	sum[gmax], sumsum, sums[gmax], sumsums, sumu[gmax], sumsumu;
	
	printf("Counting rooted");
	if (rt==1) printf(" and sensed");
	if (rt>=2) printf(", sensed and unsensed");
	printf(" maps.\n");
	while (1)
	{
		printf("\nEnter maximum number of edges <= %d or -1 to quit: ",emax);
		scanf("%d",&emx);
		if (emx<0) break;
		if (emx>emax) continue;
		printf("Enter minimum number of edges between 0 and %d: ",emx);
		scanf("%d",&emn);
		if (emn>emx) emn=emx;
		if (emn<0) emn=0;
		for (e=emn; e<=emx; e++)
		{
			printf("\nNumber of edges = %d\n",e);
			printf("Enter maximum number of vertices between 1 and %d: ",e+1);
			scanf("%d",&vmx);
			if (vmx>e+1) vmx=e+1;
			printf("Enter minimum number of vertices between 1 and %d: ",vmx);
			scanf("%d",&vmn);
			printf("Enter 1 to print the Lehman code for each rooted map or 0 not to: ");
			scanf("%d",&prt);
			printf("\nNumber of rooted maps with e edges and v vertices.\n");
			printf("  e   v");
			for (g=0;g<=e/2;g++) printf("         g=%1d",g);
			printf("       all g\n");
			sumsum=0;
			for (g=0;g<=e/2;g++) sum[g]=0;
			if (rt>0)
			{
				sumsums=0;
				for (g=0;g<=e/2;g++) sums[g]=0;
			}
			if (rt>1)
			{
				sumsumu=0;
				for (g=0;g<=e/2;g++) sumu[g]=0;
			}
			for (v=vmn; v<=vmx; v++)
			{
				np=v-1;
				nq=e-np;
				for (g=0;g<=e/2;g++) root[g]=0;
				ct=0;
				if (rt>0)
				{
					cts[v]=0;
					for (g=0;g<=e/2;g++) sensed[v][g]=0;
				}
				if (rt>1)
				{
					ctu[v]=0;
					for (g=0;g<=e/2;g++) unsensed[v][g]=0;
				}
				firstint(nq);
				do
				{
					qpairs(nq);
					g=genus(nq);
					firstpar(np);
					do
					{
						if (ppairs(np)>0) {printf("Invalid parenthesis system.\n"); break;}
						firstsys(np,nq);
						do
						{
							ct++;
							root[g]++;
							sumsum++;
							sum[g]++;
							if (prt)
							{
								printf("%5ld:",ct);
								for (i=1; i<=2*e; i++)
								{
									if (s[i]==0) printf("  (");
									else if (s[i]==-1) printf("  )");
									else printf("%3d",s[i]);
								}
								printf("  genus = %d",g);
							}
							if (rt>0)
							{
								decode(e,v,vert,next);
								max=1;
								for (r=2;r<=2*e;r++)
								{
									i = Lehman(e, v, vert, next, r, code, 1);
									if (i>0)
									{
										max=0;
										break;
									}
								}
								if (max)
								{
									cts[v]++;
									sensed[v][g]++;
									sumsums++;
									sums[g]++;
									if (prt) printf(" Sensed map #%ld.",sensed[v][g]);
									if (rt>1)
									{
										for (i=1;i<=2*e;i++) prev[next[i]]=i;
										for (r=1;r<=2*e;r++)
										{
											i = Lehman(e, v, vert, prev, r, code, 1);
											if (i>0)
											{
												max=0;
												break;
											}
										}
										if (max)
										{
											ctu[v]++;
											unsensed[v][g]++;
											sumsumu++;
											sumu[g]++;
											if (prt) printf(" Unsensed map #%ld",unsensed[v][g]);
										} /* end if reflex-maximality test */
									} /* end if rt>1 */
								} /* end if maximality test */
							} /* end if rt>0 */
							if (prt) printf("\n");
						}
						while (nextsys(np,nq)>0);
					}
					while (nextpar(np)>0);
				}
				while (nextint(nq)>0);
				printf("%3d%4d",e,v);
				for (g=0;g<=e/2;g++) printf("%12ld",root[g]);
				printf("%12ld\n",ct);
			} /* end for v */
			printf("%3d sum",e);
			for (g=0;g<=e/2;g++) printf("%12ld",sum[g]);
			printf("%12ld\n",sumsum);
			if (rt>0)
			{
				printf("\nNumber of sensed maps with e edges and v vertices.\n");
				printf("  e   v");
				for (g=0;g<=e/2;g++) printf("         g=%1d",g);
				printf("       all g\n");
				for (v=vmn; v<=vmx; v++)
				{
					printf("%3d%4d",e,v);
					for (g=0;g<=e/2;g++) printf("%12ld",sensed[v][g]);
					printf("%12ld\n",cts[v]);
				}
				printf("%3d sum",e);
				for (g=0;g<=e/2;g++) printf("%12ld",sums[g]);
				printf("%12ld\n",sumsums);				
			} /* end if rt>0 */
			if (rt>1)
			{
				printf("\nNumber of unsensed maps with e edges and v vertices.\n");
				printf("  e   v");
				for (g=0;g<=e/2;g++) printf("         g=%1d",g);
				printf("       all g\n");
				for (v=vmn; v<=vmx; v++)
				{
					printf("%3d%4d",e,v);
					for (g=0;g<=e/2;g++) printf("%12ld",unsensed[v][g]);
					printf("%12ld\n",ctu[v]);
				}
				printf("%3d sum",e);
				for (g=0;g<=e/2;g++) printf("%12ld",sumu[g]);
				printf("%12ld\n",sumsumu);				
			} /* end if rt>1 */			
		} /* end for e*/
	} /* end while */
	printf("End of execution of count.\n");
} /* end count */

void pcount(int rt)
{
	int i, v, e, np, nq, emx, emn, vmx, vmn, prt, r, vert[dim], next[dim], code[dim], prev[dim], max;
	long int ct, root, sensed[dim], unsensed[dim], sum, sums, sumu;
	
	printf("Counting rooted");
	if (rt==1) printf(" and sensed");
	if (rt>=2) printf(", sensed and unsensed");
	printf(" planar maps.\n");
	while (1)
	{
		printf("\nEnter maximum number of edges <= %d or -1 to quit: ",emax);
		scanf("%d",&emx);
		if (emx<0) break;
		if (emx>emax) continue;
		printf("Enter minimum number of edges between 0 and %d: ",emx);
		scanf("%d",&emn);
		if (emn>emx) emn=emx;
		if (emn<0) emn=0;
		for (e=emn; e<=emx; e++)
		{
			printf("\nNumber of edges = %d\n",e);
			printf("Enter maximum number of vertices between 1 and %d: ",e+1);
			scanf("%d",&vmx);
			if (vmx>e+1) vmx=e+1;
			printf("Enter minimum number of vertices between 1 and %d: ",vmx);
			scanf("%d",&vmn);
			printf("Enter 1 to print the Lehman code for each rooted map or 0 not to: ");
			scanf("%d",&prt);
			printf("\nNumber of rooted planar maps with e edges and v vertices.\n");
			printf("  e   v\n");
			sum=0;
			if (rt>0) sums=0;
			if (rt>1) sumu=0;
			for (v=vmn; v<=vmx; v++)
			{
				np=v-1;
				nq=e-np;
				root=0;
				ct=0;
				if (rt>0) sensed[v]=0;
				if (rt>1) unsensed[v]=0;
				firstbrac(nq);
				do
				{
					if (bpairs(nq)>0) {printf("Invalid bracket system.\n"); break;}
					firstpar(np);
					do
					{
						if (ppairs(np)>0) {printf("Invalid parenthesis system.\n"); break;}
						firstsys(np,nq);
						do
						{
							ct++;
							root++;
							sum++;
							if (prt)
							{
								printf("%5ld:",ct);
								for (i=1; i<=2*e; i++)
								{
									if (s[i]==0) printf(" (");
									else if (s[i]==-1) printf(" )");
									else if (s[i]==2) printf(" [");
									else if (s[i]==1) printf(" ]");
									else {printf(" Oops!\n"); return;}
								}
							}
							if (rt>0)
							{
								decode(e,v,vert,next);
								max=1;
								for (r=2;r<=2*e;r++)
								{
									i = Alfred(e, v, vert, next, r, code, 1);
									if (i>0)
									{
										max=0;
										break;
									}
								}
								if (max)
								{
									sensed[v]++;
									sums++;
									if (prt) printf("  Sensed map #%ld.",sensed[v]);
									if (rt>1)
									{
										for (i=1;i<=2*e;i++) prev[next[i]]=i;
										for (r=1;r<=2*e;r++)
										{
											i = Alfred(e, v, vert, prev, r, code, 1);
											if (i>0)
											{
												max=0;
												break;
											}
										}
										if (max)
										{
											unsensed[v]++;
											sumu++;
											if (prt) printf(" Unsensed map #%ld.",unsensed[v]);
										} /* end if reflex-maximality test */
									} /* end if rt>1 */
								} /* end if maximality test */
							} /* end if rt>0 */
							if (prt) printf("\n");
						}
						while (nextsys(np,nq)>0);
					}
					while (nextpar(np)>0);
				}
				while (nextbrac(nq)>0);
				printf("%3d%4d",e,v);
				printf("%12ld\n",root);
			} /* end for v */
			printf("%3d sum",e);
			printf("%12ld\n",sum);
			if (rt>0)
			{
				printf("\nNumber of sensed planar maps with e edges and v vertices.\n");
				printf("  e   v\n");
				for (v=vmn; v<=vmx; v++)
				{
					printf("%3d%4d",e,v);
					printf("%12ld\n",sensed[v]);
				}
				printf("%3d sum",e);
				printf("%12ld\n",sums);				
			} /* end if rt>0 */
			if (rt>1)
			{
				printf("\nNumber of unsensed planar maps with e edges and v vertices.\n");
				printf("  e   v\n");
				for (v=vmn; v<=vmx; v++)
				{
					printf("%3d%4d",e,v);
					printf("%12ld\n",unsensed[v]);
				}
				printf("%3d sum",e);
				printf("%12ld\n",sumu);				
			} /* end if rt>1 */			
		} /* end for e*/
	} /* end while */
	printf("End of execution of pcount.\n");
} /* end pcount */

void test(int rt)
{
	int i,d,np,nq,flag,r,max,occ[emaxp1],vert[dim],next[dim],code[dim],prev[dim];
	
	if (rt==0) printf("Testing canonical mixing.\n");
	else printf("Testing Lehman coding and decoding and maximality.\n");
	if (rt>1) printf("Also testing maximality compared to reflected map.\n");
	while (1) 
	{
		printf("Enter number of integer pairs <= %d, -1 to quit: ",emax);
		scanf("%d",&nq);
		if (nq<0) break;
		if (nq>emax) continue;
		printf("Enter each integer from 1 to %d twice with the first occurrences in increasing order.\n",nq);
		flag=0;
		for (i=1;i<=nq;i++) occ[i]=0;
		for (i=0;i<2*nq;i++)
		{
			scanf("%d",&q[i]);
			occ[q[i]]++;
			if (occ[q[i]]==1) qpos1[q[i]]=i; else if (occ[q[i]]==2) qpos2[q[i]]=i; else flag=1;
		}
		q[2*nq]=0;
		for (i=1; i<=nq; i++)
			if (occ[i]!=2) {flag=1; break;}
		if (qpos1[1]!=0)  flag=1;
		for (i=2; i<=nq; i++)
			if (qpos1[i]<=qpos1[i-1]) {flag=1; break;}
		if (flag) {printf("Invalid integer system\n"); continue;}
		qpairs(nq);
		for (i=0; i<2*nq; i++) printf("%3d",q[i]);
		printf("   ");
		for (i=0; i<2*nq; i++) printf("%3d",qmate[i]);
		printf("  genus = %d\n", genus(nq));
		printf("Enter number of parenthesis pairs <=%d, -1 to quit: ",emax-nq);
		scanf("%d",&np);
		printf("Enter balanced parenthesis system with %d pairs, (=0, )=-1.\n",np);
		for (i=1;i<=2*np;i++)
		{
			scanf("%d",&p[i]);
			if ((p[i]!=0)&&(p[i]!=-1)) flag=1;
		}
		p[0]=-1;
		if (!flag) flag=ppairs(np);
		if (flag) {printf("Invalid parenthesis system.\n"); continue;}
		for (i=1; i<=2*np; i++) printf("%3d",p[i]);
		printf("   ");
		for (i=1; i<=2*np; i++) printf("%3d",pmate[i]);
		printf("\n");
		printf("Enter 1 to print details of the calculations, 0 not to: ");
		scanf("%d",&bug);
		firstsys(np,nq);
		do
		{
			for (i=1; i<=2*(np+nq); i++) printf("%3d",s[i]);
			printf("   ");
			for (i=1; i<=2*(np+nq); i++) printf("%3d",smate[i]);
			printf("\n");
			if (rt>0)
			{
				printf("Enter 1 to decode this system or 0 not to: ");
				scanf("%d",&flag);
				if (flag)
				{
					decode(np+nq, np+1, vert, next);
					printf("vert: ");
					for (i=1; i<=2*(np+nq); i++) printf("%3d",vert[i]);
					printf("\n");
					printf("next: ");
					for (i=1; i<=2*(np+nq); i++) printf("%3d",next[i]);
					printf("\n");
					i = Lehman(np+nq, np+1, vert, next, 1, code, 0);
					printf("And here is the code of the rooted map.\n");
					for (i=1;i<=2*(np+nq);i++)
					{
						printf("%3d",code[i]);
						if (code[i]!=s[i])
						{
							printf("  Oops!\n"); break;
						}
					}
					printf("\n");
					printf("Enter 1 to test this code for maximality or 0 not to: ");
					scanf("%d",&flag);
					if (flag)
					{
						max=1;
						for (r=2; r<=2*(np+nq); r++)
						{
							d = Lehman(np+nq, np+1, vert, next, r, code, 1);
							if (d>0)
							{
								printf("This code is not maximal.  Rooting at %d is bigger.\n",r);
								i = Lehman(np+nq, np+1, vert, next, r, code, 0);
								for (i=1; i<=2*(np+nq); i++) printf("%3d", code[i]);
								printf("\n");
								max=0;
								break;
							}
						}
						if (max)
						{
							printf("This code is maximal.\n");
							if (rt>1)
							{
								for (i=1;i<=2*(np+nq);i++) prev[next[i]]=i;
								printf("vert: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",vert[i]);
								printf("\n");
								printf("next: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",next[i]);
								printf("\n");
								printf("prev: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",prev[i]);
								printf("\n");
								for (r=1; r<=2*(np+nq); r++)
								{
									d = Lehman(np+nq, np+1, vert, prev, r, code, 1);
									if (d>0)
									{
										printf("This code is not reflect-maximal.  Rooting at %d is bigger.\n",r);
										i = Lehman(np+nq, np+1, vert, prev, r, code, 0);
										for (i=1; i<=2*(np+nq); i++) printf("%3d", code[i]);
										printf("\n");
										max=0;
										break;
									}
								}
								if (max) printf("This code is reflect-maximal.\n");
							} /* end if reflect-maximality test flag */
						} /* end if max */
					} /* end if maximality test flag */
				} /* end if coding flag */
			} /* end if rt>0 */
		}
		while (nextsys(np,nq)>0);
	}
	printf("End of execution of test.\n");
}

void ptest(int rt)
{
	int i,d,np,nq,flag,r,max,vert[dim],next[dim],code[dim],prev[dim];
	
	if (rt==0) printf("Testing canonical mixing of pairs of Dyck words.\n");
	else printf("Testing Lehman coding and decoding and maximality.\n");
	if (rt>1) printf("Also testing maximality compared to reflected map.\n");
	while (1) 
	{
		printf("Enter number of bracket pairs <= %d, -1 to quit: ",emax);
		scanf("%d",&nq);
		if (nq<0) break;
		if (nq>emax) continue;
		printf("Enter balanced bracket system with %d pairs, [=2, ]=1.\n",nq);
		flag=0;
		for (i=0;i<2*nq;i++)
		{
			scanf("%d",&q[i]);
			if ((q[i]!=2)&&(q[i]!=1)) flag=1;
		}
		q[2*nq]=0;
		if (!flag) flag=bpairs(nq);
		if (flag) {printf("Invalid bracket system.\n"); continue;}
		for (i=0; i<2*nq; i++) printf("%3d",q[i]);
		printf("   ");
		for (i=0; i<2*nq; i++) printf("%3d",qmate[i]);
		printf("\n");
		printf("Enter number of parenthesis pairs <=%d, -1 to quit: ",emax-nq);
		scanf("%d",&np);
		printf("Enter balanced parenthesis system with %d pairs, (=0, )=-1.\n",np);
		for (i=1;i<=2*np;i++)
		{
			scanf("%d",&p[i]);
			if ((p[i]!=0)&&(p[i]!=-1)) flag=1;
		}
		p[0]=-1;
		if (!flag) flag=ppairs(np);
		if (flag) {printf("Invalid parenthesis system.\n"); continue;}
		for (i=1; i<=2*np; i++) printf("%3d",p[i]);
		printf("   ");
		for (i=1; i<=2*np; i++) printf("%3d",pmate[i]);
		printf("\n");
		printf("Enter 1 to print details of the calculations, 0 not to: ");
		scanf("%d",&bug);
		firstsys(np,nq);
		do
		{
			for (i=1; i<=2*(np+nq); i++) printf("%3d",s[i]);
			printf("   ");
			for (i=1; i<=2*(np+nq); i++) printf("%3d",smate[i]);
			printf("\n");
			if (rt>0)
			{
				printf("Enter 1 to decode this system or 0 not to: ");
				scanf("%d",&flag);
				if (flag)
				{
					decode(np+nq, np+1, vert, next);
					printf("vert: ");
					for (i=1; i<=2*(np+nq); i++) printf("%3d",vert[i]);
					printf("\n");
					printf("next: ");
					for (i=1; i<=2*(np+nq); i++) printf("%3d",next[i]);
					printf("\n");
					i = Alfred(np+nq, np+1, vert, next, 1, code, 0);
					printf("And here is the code of the rooted map.\n");
					for (i=1;i<=2*(np+nq);i++)
					{
						printf("%3d",code[i]);
						if (code[i]!=s[i])
						{
							printf("  Oops!\n"); break;
						}
					}
					printf("\n");
					printf("Enter 1 to test this code for maximality or 0 not to: ");
					scanf("%d",&flag);
					if (flag)
					{
						max=1;
						for (r=2; r<=2*(np+nq); r++)
						{
							d = Alfred(np+nq, np+1, vert, next, r, code, 1);
							if (d>0)
							{
								printf("This code is not maximal.  Rooting at %d is bigger.\n",r);
								i = Alfred(np+nq, np+1, vert, next, r, code, 0);
								for (i=1; i<=2*(np+nq); i++) printf("%3d", code[i]);
								printf("\n");
								max=0;
								break;
							}
						}
						if (max)
						{
							printf("This code is maximal.\n");
							if (rt>1)
							{
								for (i=1;i<=2*(np+nq);i++) prev[next[i]]=i;
								printf("vert: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",vert[i]);
								printf("\n");
								printf("next: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",next[i]);
								printf("\n");
								printf("prev: ");
								for (i=1; i<=2*(np+nq); i++) printf("%3d",prev[i]);
								printf("\n");
								for (r=1; r<=2*(np+nq); r++)
								{
									d = Alfred(np+nq, np+1, vert, prev, r, code, 1);
									if (d>0)
									{
										printf("This code is not reflect-maximal.  Rooting at %d is bigger.\n",r);
										i = Alfred(np+nq, np+1, vert, prev, r, code, 0);
										for (i=1; i<=2*(np+nq); i++) printf("%3d", code[i]);
										printf("\n");
										max=0;
										break;
									}
								}
								if (max) printf("This code is reflect-maximal.\n");
							} /* end if reflect-maximality test flag */
						} /* end if max */
					} /* end if maximality test flag */
				} /* end if coding flag */
			} /* end if rt>0 */
		}
		while (nextsys(np,nq)>0);
	} /* end while(1) */
	printf("End of execution of test.\n");
}

int Lehman(int e, int v, int vert[dim], int next[dim], int r, int code[dim], int c)
/* Constructs the Lehman code for the map defined by vert and next rooted at r. */
{
	int d, i, el[dim], vl[dim], ec;
	
	for (i=1; i<=2*e; i++) el[i]=0;
	for (i=1; i<=v; i++) vl[i]=0;
	d=r; ec=0; vl[vert[d]]=1;
	for (i=1; i<=2*e; i++)
	{
		if (el[d]==0) /* The edge containing d is seen for the first time. */
		{
			if (vl[vert[2*e+1-d]]==0) /* The vertex at the other end of this edge not yet seen. */
			{
				el[d]=-1; el[2*e+1-d]=-1; /* Put this edge in the tree. */
				code[i]=0; /* left parenthesis */
				vl[vert[2*e+1-d]]=1; /* Now we've seen this vertex. */
			}
			else /* The vertex at the other end of this edge has already been seen. */
			{
				ec++;
				el[d]=ec; el[2*e+1-d]=ec; /*Put this edge in the cotree.*/
				code[i]=ec;
			}
		}
		else code[i]=el[d]; /* The edge containing d is seen for the second time. */
		if (c>0) /* We're testing the code s, decoded by vert and next, for maximality. */
		{
			if (code[i]>s[i]) return d; /* The map rooted at d has a bigger code than s. */
			else if (code[i]<s[i]) return 0; /* The map rooted at d has a smaller code than s. */
		}
		if (code[i]<=0) d=next[2*e+1-d]; /* Go to the other end of the edge containing d and rotate. */
		else d=next[d]; /* Just rotate. */
	} /* end for */
	return 0;
}

int Alfred(int e, int v, int vert[dim], int next[dim], int r, int code[dim], int c)
/* Constructs the Lehman code for the planar map defined by vert and next rooted at r. */
{
	int d, i, el[dim], vl[dim];
	
	for (i=1; i<=2*e; i++) el[i]=0;
	for (i=1; i<=v; i++) vl[i]=0;
	d=r; vl[vert[d]]=1;
	for (i=1; i<=2*e; i++)
	{
		if (el[d]==0) /* The edge containing d is seen for the first time. */
		{
			if (vl[vert[2*e+1-d]]==0) /* The vertex at the other end of this edge not yet seen. */
			{
				el[d]=-1; el[2*e+1-d]=-1; /* Put this edge in the tree. */
				code[i]=0; /* left parenthesis */
				vl[vert[2*e+1-d]]=1; /* Now we've seen this vertex. */
			}
			else /* The vertex at the other end of this edge has already been seen. */
			{
				el[d]=1; el[2*e+1-d]=1; /*Put this edge in the cotree.*/
				code[i]=2;
			}
		}
		else code[i]=el[d]; /* The edge containing d is seen for the second time. */
		if (c>0) /* We're testing the code s, decoded by vert and next, for maximality. */
		{
			if (code[i]>s[i]) return d; /* The map rooted at d has a bigger code than s. */
			else if (code[i]<s[i]) return 0; /* The map rooted at d has a smaller code than s. */
		}
		if (code[i]<=0) d=next[2*e+1-d]; /* Go to the other end of the edge containing d and rotate. */
		else d=next[d]; /* Just rotate. */
	} /* end for */
	return 0;
}

void decode(int e, int v, int vert[dim], int next[dim])
/* Creates the map defined by vert and next rooted at dart 1 whose Lehman code is s. */
{
	int vc, ec, dold, d, i, ds[dim];
	
	vc=1; ec=0; dold=0; vert[0]=1;
	for (i=1; i<=2*e; i++) ds[i]=0;
	for (i=1; i<=2*e; i++)
	{
		if (ds[i]==0) /* First occurrence of an integer or a left parenthesis. */
		{
			ec++; d=ec; ds[i]=d; ds[smate[i]]=2*e+1-d;
		}
		else d=ds[i];
		vert[d]=vert[dold];
		next[dold]=d;
		if (s[i]<=0) dold=2*e+1-d;
		else dold=d;
		if (s[i]==0) { vc++; vert[dold]=vc;}
	}
	vert[1]=1; next[dold]=1;
}

int nextsys(int np, int nq)
{
	int ipr, ipos, jpos, kpos, t, i, ex, sti[emax], rt[dim], lf[dim];
	
	for (i=0; i<=2*(np+nq); i++) { rt[i]=i+1; lf[i]=i-1; }
	rt[2*(np+nq)]=0; lf[0]=2*(np+nq);
	ipr=2*np;
	while (ipr>0)
	{
		ipos=prpos[ipr]; /* the position of the parentheses having priority ipr */
		jpos=rt[ipos];
		t=s[jpos];
		if ((jpos!=0)&&(t>0)&&((s[ipos]==0)||(smate[ipos]<=smate[jpos]))) /* the parenthesis can be moved to the right */
		{
			if (bug) printf("The parenthesis %d in position %d can be moved.\n",s[ipos],ipos);
			s[ipos+1]=s[ipos];
			s[ipos]=t;
			if (jpos>ipos+1)
			{
				smate[ipos+1]=smate[ipos];
				smate[smate[ipos]]++;
				smate[ipos]=smate[jpos];
				smate[smate[jpos]]=ipos;					
				s[jpos]=-2;
				rt[ipos]=ipos+1;
				lf[ipos+1]=ipos;
				rt[ipos+1]=rt[ipos];
				lf[rt[jpos]]=ipos+1;
			}
			else
			{
				ex=smate[ipos];
				smate[ipos]=smate[jpos];
				smate[jpos]=ex;
				smate[smate[ipos]]=ipos;
				smate[smate[jpos]]=jpos;
			}
			prpos[ipr]++;
			if (bug)
			{	
				for (i=1; i<=2*(np+nq); i++) printf("%3d",s[i]);
				printf("   ");
				for (i=1; i<=2*(np+nq); i++) printf("%3d",smate[i]);
				printf("\n");
			}
			/* now move the integers to replace the deleted parentheses */
			jpos=2*(np+nq); kpos=jpos-1;
			while (kpos>=1)
			{
				if (kpos==jpos) kpos--;
				else if ((s[kpos]==0)||(s[kpos]==-1)) {jpos=kpos; kpos--;}
				else if (s[kpos]==-2)
					{
						if (s[jpos]>1) jpos--;
						kpos--;
					}
				else /* kpos<jpos and s[kpos] is an integer */
				{
					if (s[jpos]==-2) /* it's a blank - insert the integer */
					{
						s[jpos]=s[kpos]; s[kpos]=-2;
						smate[jpos]=smate[kpos];
						smate[smate[jpos]]=jpos;
						kpos--;
					} /* end of inserting the integer */
					jpos--;
				} /* end else - kpos<jpos and s[kpos] is an integer */
			} /* end while - all blanks are now contiguous with nearest parenthesis to their left */
			if (bug)
			{
				printf("Now the integers have been moved into the blanks.\n");
				for (i=1; i<=2*(np+nq); i++) printf("%3d",s[i]);
				printf("   ");
				for (i=1; i<=2*(np+nq); i++) printf("%3d",smate[i]);
				printf("\n");
			}
			/* now insert the deleted parenthesis in the blanks */
			i=0;
			for (ipos=1; ipos<=2*(np+nq); ipos++)
			{
				if (s[ipos]<=0)
				{
					i++;
					if (s[ipos]==-2)
					{
						s[ipos]=p[i];
						prpos[pr[i]]=ipos;
					}
				}
			} /* end for - now the deleted parentheses have been inserted */
			/* now find the mates to all the parentheses */
			ex=-1;
			for (i=1; i<=2*(np+nq); i++)
			{
				if (s[i]==0) /* it's a left parenthesis - save its position */
				{ ex++; sti[ex]=i; }
				else if (s[i]==-1) /* it's a right parenthesis */
				{ /* the position of its mate is on the top of the stack */
					smate[i]=sti[ex];
					smate[sti[ex]]=i;
					ex--;
				}
			}
			if (bug)
			{
				printf("And now the deleted parentheses have been replaced.\n");
				for (i=1; i<=2*(np+nq); i++) printf("%3d",s[i]);
				printf("   ");
				for (i=1; i<=2*(np+nq); i++) printf("%3d",smate[i]);
				printf("\n");
			}
			return 1;
		} /* end if the parenthesis can be moved to the right */
		else
		{
			if (bug) printf("The parenthesis %d in position %d will be deleted.\n",s[ipos],ipos);
			s[ipos]=-2; /* it can't be moved - delete it */
			if (bug)
			{
				for (i=1; i<=2*(np+nq); i++) printf("%3d",s[i]);
				printf("   ");
				for (i=1; i<=2*(np+nq); i++) printf("%3d",smate[i]);
				printf("\n");
			}
			rt[lf[ipos]]=rt[ipos];
			lf[rt[ipos]]=lf[ipos];
			ipr--;
		}
	} /* end while ipr>1 - it was already the last system */
	return 0;
} /* end of nextsys */

void firstsys(int np, int nq)
{
	int i,j;
	
	i=0;
	while (i<2*np)
	{
		i++;
		s[i]=p[i];
		smate[i]=pmate[i];
	}
	j=0;
	while (j<2*nq)
	{
		i++;
		s[i]=q[j];
		smate[i]=qmate[j]+2*np+1;
		j++;
	}
}

int ppairs(int n)
{
	int i, ex, prc, sti[emax], stpr[emax];
	
	ex=-1;
	prc=-1;
	for (i=1; i<=2*n; i++)
	{
		if (p[i]==0)
		{
			ex++;
			sti[ex]=i;
			prc++; prc++;
			prpos[prc]=i;
			pr[i]=prc;
			stpr[ex]=prc;
		}
		else
		{
			pmate[i]=sti[ex]; pmate[sti[ex]]=i;
			prpos[stpr[ex]+1]=i;
			pr[i]=stpr[ex]+1;
			ex--;
			if (ex<-1) return 1;
		}
	}
	if (ex!=-1) return 1; else return 0;
}

void qpairs(int n)
{
	int i;
	for (i=1; i<=n; i++)
	{
		qmate[qpos1[i]]=qpos2[i];
		qmate[qpos2[i]]=qpos1[i];
	}
}

int genus(int n)
{
	int i, j, f, t, new[emax];
	
	f=0;
	for (i=0; i<2*n; i++) new[i]=1;
	for (i=0; i<2*n; i++) if (new[i])
	{
		j=i;
		do
		{
			new[j]=0;
			t=q[j];
			if (j==qpos2[t]) j=qpos1[t]; else j=qpos2[t];
			if (j==2*n-1) j=0; else j++;
		}
		while (j!=i);
		f++;
	}
	return (1+n-f)/2;
}
		
void firstint(int n)
{
	int i;
	q[2*n]=0;
	for (i=1; i<=n; i++)
	{
		q[2*i-2]=i; qpos1[i]=2*i-2;
		q[2*i-1]=i; qpos2[i]=2*i-1;
		dir[i]=1;
	}
}

int nextint(int n)
{
	int i,t,qp;
	
	i=1;
	while(i<n)
	{
		qp=qpos2[i];
		if (dir[i]>0)
		{
			t=q[qp+1];
			if (t<i) {dir[i]=-1; i++;}
			else
			{
				q[qp+1]=i;
				q[qp]=t;
				qpos2[i]++;
				if (qpos2[i]==qpos1[t]) qpos1[t]--; else qpos2[t]--;
				return 1;
			}
		}
		else
		{
			t=q[qp-1];
			if (t==i) {dir[i]=1; i++;}
			else
			{
				q[qp-1]=i;
				q[qp]=t;
				qpos2[i]--;
				if (qpos2[i]==qpos1[t]) qpos1[t]++; else qpos2[t]++;
				return 1;
			}
		}
	}
	return 0;
}

void firstpar(int n)
{
	int i;
	
	p[0]=-1;
	i=0;
	while (i<2*n) {i++; p[i]=0; i++; p[i]=-1;}
}

int nextpar(int n)
{
	int i,exc;

	if (n<=0) return 0;
	i=2*n;
    exc=0;
	while (p[i]==-1) {i--; exc++;}
	while (p[i]==0) {i--; exc--;}
	if (i==0) return 0;
	p[i]=0;
	exc++; exc++;
	while (exc>0) {i++; p[i]=-1; exc--;}
	while (i<2*n) {i++; p[i]=0; i++; p[i]=-1;}
	return 1;
}

void firstbrac(int n)
{
	int i;
	
	i=-1;
	while (i<2*n-1) {i++; q[i]=2; i++; q[i]=1;}
}

int nextbrac(int n)
{
	int i,exc;
	
	if (n<=0) return 0;
	i=2*n-1;
    exc=0;
	while (q[i]==1) {i--; exc++;}
	while ((i>=0)&&(q[i]==2)) {i--; exc--;}
	if (i<0) return 0;
	q[i]=2;
	exc++; exc++;
	while (exc>0) {i++; q[i]=1; exc--;}
	while (i<2*n-1) {i++; q[i]=2; i++; q[i]=1;}
	return 1;
}

int bpairs(int n)
{
	int i, ex, prc, sti[emax], stpr[emax];
	
	ex=-1;
	prc=-1;
	for (i=0; i<2*n; i++)
	{
		if (q[i]==2)
		{
			ex++;
			sti[ex]=i;
			prc++; prc++;
			stpr[ex]=prc;
		}
		else
		{
			qmate[i]=sti[ex]; qmate[sti[ex]]=i;
			pr[i]=stpr[ex]+1;
			ex--;
			if (ex<-1) return 1;
		}
	}
	if (ex!=-1) return 1; else return 0;
}